Si vuole calcolare il numero di riproduzione di base $R$ nel tempo $t$ dato il numero noto di nuovi casi $k$ in $t$.
Si userà la seguente notazione:
$X \in (a, b)$: $a < X < b$
$X \in [a, b]$: $a \leq X \leq b$
$X \in (a, b]$: $a < X \leq b$
$X \in [a, b)$: $a \leq X < b$
$X_i$: valore della variabile $X$ al tempo $i$
$R_t$: valore di $R$ in $t$
$k_t$: valore di $k$ in $t$
$P(X)$: probabilità a priori che si verifichi $X$
$P(X|Y^*)$: probabilità a posteriori che si verifichi $X$ noto il valore di $Y^*$
$\mathscr{L}(Y^*|X)$: verosimiglianza (likelihood) del verificarsi noto di $Y^*$ data la variabile aleatoria $X$
$\langle X |$: vettore riga $\{X_0 \cdots X_n \}$
$| X \rangle$: vettore colonna $\{X_0 \cdots X_n \}^{\mathbf{T}}$
$\big| \langle X | \big|$ o $\big| | X \rangle \big|$: modulo del vettore $X$
$\overline{X}_{n,m}$: matrice $X$ di $n$ righe e $m$ colonne
$\mathscr{N}(\mu,\sigma)$: distributione normale (Gaussiana) avente media $\mu$ e deviazione standard $\sigma$
$\mathscr{P}(X,Y)$: distribuzione di Poisson di $X$ su variabile aleatoria $Y$
Si procede calcolando la probabilità a posteriori del verificarsi di $R_t$ dato $k_t^*$ in $t$:
$$ P(R_t|k_t^*) $$dove con $x^*$ si indica che è l'evento noto (già accaduto).
Per il Teorema di Bayes:
$$ P(A|B^*) = \frac{P(B^*|A) \cdot P(A)}{P(B^*)} $$quindi
$$ P(R_t|k_t^*) = \frac{P(k_t^*|R) \cdot P(R_t)}{P(k_t^*)} $$sapendo che
$$ P(k_t^*) = \sum \bigg( P(k_t|R_t) \cdot P(R_t) \bigg) $$calcolando dunque, per ogni $t$:
$P(k_t^*|R)$: la probabilità a posteriori di aver osservato $k_t^*$ nuovi casi dati i possibili valori di $R$, moltiplicato
$P(R_t)$: la probabilità a priori dei possibili valori di $R_t$ in $t$ in base al valore calcolato del giorno precedente $R_{t-1}$, diviso
$P(k_t^*)$: la probabilità a priori di osservare $k_t^*$ nuovi casi in generale (pari alla sommatoria dei numeratori precedenti)
È noto che $R_t$ dipende da $\gamma$, l'inverso dell'intervallo seriale $T_{serial}$ ovvero il tempo medio tra casi successivi in una catena di trasmissione, spesso calcolato come intervallo tra presunta infezione e comparsa dei sintomi (periodo di contagiosità asintomatico). Per Covid-19 è stato stimato intorno ai $7.5$ giorni:
Si crea un vettore riga di tutti i possibili valori di $R$ da $R_{min}$, pari a $0$, ad $R_{max}$ per il quale si può assumere il valore massimo $12$.
Sarà utile anche un vettore colonna dei possibili valori di $R$.
$$ | R \rangle = \langle R |^{\mathbf{T}} $$$$ | R \rangle = \begin{Bmatrix} R_{min}\\ \vdots\\ R_{max} \end{Bmatrix}_{1201,1} $$Si crea un vettore riga $\langle k_t |$ di tutti i $k$ nuovi casi per giorno da $t=0$ fino a $t=\omega$ (oggi).
Per risolvere $P(k_t^*|R)$ si nota anzitutto che, data una $\lambda$ serie di possibili $k$ nuovi casi
$$ P(k|\lambda) = \mathscr{P}(k, \lambda) = \frac{\lambda^k \cdot e^{-\lambda}}{k!}$$ovvero, la probabilità a posteriori di avere $k$ nuovi casi dati $\lambda$ possibili valori $k$ è una distribuzione di Poisson.
Sapendo però che il valore di $k$ è noto (tutti i valori di $k^*$ sono nel vettore riga $\langle k_t |$) si può dire che
$$ P(k^*|\lambda) \propto \mathscr{L}(\lambda|k^*) $$ovvero, la probabilità di aver avuto $k^*$ nuovi casi date $\lambda$ possibilità è proporzionale alla verosimiglianza delle $\lambda$ possibilità noto $k^*$.
Dalla letteratura, è nota la relazione tra $\lambda$, $k_t$ ed $R_t$:
$$ \lambda = k_t \cdot \exp{\big[ \gamma (R_t - 1) \big]} $$Grazie a questa relazione è possibile dire che
$$ P(k_t^*|R) = P(k_t^*|\lambda) $$e costruire una matrice $\overline{\lambda}$ di $\omega-1$ righe (perché il primo giorno è di outcome) e $1201$ colonne contenente per ogni giorno $t$ tutti i possibili valori di $\lambda$ dati i possibili $1201$ valori di $R$.
Noti i valori della matrice $\overline{\lambda}$ si può dunque applicare la distribuzione di Poisson al vettore riga $\langle k_t |$ con $t \in (0,\omega]$ e ottenere così la matrice $\overline{P(k_t^*|R)}$ di tutte le verosimiglianze dei possibili $1201$ valori di $R$ dati gli $\omega-1$ valori di $k_t^*$:
La probabilità a priori dei possibili valori $R$ può essere interpretata come una ditribuzione normale (Gaussiana) avente media $\mu = R$ e deviazione standard $\sigma$ che può stimata a $0.25$.
Si può dunque creare una matrice $\overline{P(R)}$ di $1201 \times 1201$ i cui elementi $P(n)_m$ di riga $n$ e colonna $m$ sono la probabilità di $R=n$ con $\mu=m$
È dunque possibile ora creare una matrice vuota $\overline{P(R_t|k_t^*)}$ di $1201 \times \omega$ i cui elementi $P(R_n)_m$ di riga $n$ e colonna $m$ saranno i valori della probabilità di $R=m$ al tempo $t=n$.
Non avendo alcuna informazione precedente a $t=0$, si imposta la probabilità iniziale della prima colonna $| P(R_0|k_0^*) \rangle$ ovvero degli elementi $P(R_0)_m$ come costante in modo che $\sum P(R_0)_m = 1$:
Per ogni giorno $t \in (0, \omega]$ si calcola la probabilità di $R_t$ basata sulla distribuzione normale del giorno precedente $t-1$
e si moltiplica per $| P(k_t^*|R) \rangle$ calcolato in precedenza, ottenendo così il numeratore della funzione 1:
$$ | P(k_t^*|R) \rangle \cdot | P(R_t)_m \rangle $$Il denomiatore $P(k_t^*)$ è la sommatoria dei numeratori precedenti:
ottenendo dunque per ogni giorno $t$ la probabilitÃ
$$ | P(R_t|k_t^*) \rangle = \frac{ | P(k_t^*|R) \rangle \cdot | P(R_t)_m \rangle }{ \sum_{i=0}^{t} \bigg( | P(k_i^*|R) \rangle \cdot | P(R_i)_m \rangle \bigg) } $$